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Abstract 


Reconnaissance  missions  and  satellites  collect  hundreds  of  images  loaded  with 
valuable  information  to  be  utilized  by  the  Air  Force.  Intelligence  operations  must 
analyze  these  images  to  extract  the  information  needed  to  help  commanders  make 
important  decisions.  No  matter  how  obtained,  images  of  this  type  are  often  degraded  by 
noise  due  to  disruptions  such  as  atmospheric  disturbances,  optical  system  variations, 
motion,  and  large  distance  from  the  sensor  to  the  source.  This  noise  must  be  removed 
effectively  to  improve  the  quality  of  these  images  and  ensure  that  the  information 
contained  in  them  can  be  correctly  extracted.  The  Air  Force  relies  on  the  critical 
information  contained  in  these  images  and  therefore  effective  noise  removal  is  critical. 

Limitations  in  current  image  restoration  techniques  cause  visually  disturbing 
effects  in  reconstructed  images.  Therefore,  a  more  effective  algorithm  for  noise  removal 
is  required.  This  thesis  will  explore  the  use  of  redundant  wavelet-based  multiresolution 
image  processing  and  a  priori  information  to  enhance  current  image  restoration 
techniques.  Often  on  an  airplane,  tank  or  other  military  vehicle  or  building  there  may  be 
a  known  symbol,  which  can  be  used  as  perfect  a  priori  information  to  restore  the  entire 
image.  An  adaptive  redundant  discrete  wavelet  transform  will  be  trained  using  this  a 
priori  information.  Constant  and  adaptive  thresholds  are  used  on  the  noisy  image  while 
in  the  wavelet  domain.  The  quality  of  the  reconstructed  image  is  improved  over  non- 
adaptive  techniques. 


XI 


REDUNDANT  WAVELET-BASED  IMAGE  RESTORATION 


USING  A  PRIORI  INFORMATION 


I.  Introduction 


1.1  Background 

Reconnaissance  missions  and  satellites  collect  hundreds  of  images  loaded  with 
valuable  information  to  be  utilized  by  the  Air  Force.  Intelligence  operations  must 
analyze  these  images  to  extract  the  information  needed  to  help  commanders  make 
important  decisions.  No  matter  how  obtained,  images  of  this  type  are  often  degraded  by 
noise  due  to  disruptions  such  as  atmospheric  disturbances,  optical  system  variations, 
motion,  and  large  distance  from  the  sensor  to  the  source.  This  noise  must  be  removed 
effectively  to  improve  the  quality  of  these  images  and  ensure  that  the  information 
contained  in  them  can  be  correctly  extracted.  The  Air  Force  relies  on  the  critical 
information  contained  in  these  images  and  therefore  effective  noise  removal  is  critical. 


1.2  Problem 

Image  restoration  attempts  to  recover  an  image  which  has  been  degraded.  This 
thesis  will  explore  the  use  of  redundant  wavelet-based  multiresolution  image  processing 
and  a  priori  information  to  enhance  current  image  restoration  techniques.  Perfect  a  priori 
information  about  an  object  in  an  image  will  be  assumed.  Often  on  an  airplane,  tank  or 
other  military  vehicle  or  building  there  may  be  a  known  symbol,  which  can  be  used  as 
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perfect  a  priori  information  to  restore  the  entire  image.  An  adaptive  redundant  discrete 
wavelet  transform  will  be  trained  using  this  a  priori  information.  Then  the  new  adapted 
wavelet  will  be  used  to  transform  an  image  in  which  the  object  should  appear.  Wavelets 
will  be  used  because  they  have  proven  useful  in  many  signal  and  image  processing 
techniques.  Wavelet  decompositions  provide  efficient  representations  for  a  wide  range  of 
signals.  While  in  the  wavelet  domain,  the  soft  or  hard-threshold  techniques  [2]  will  be 
used  to  remove  noise.  Then  the  image  is  transformed  back  into  the  image  space  and  the 
area  of  the  image  containing  the  object  of  interest  should  be  improved.  A  successful 
algorithm  accomplishing  improved  restoration  of  an  image,  which  contains  the  training 
object,  could  be  of  great  value  to  the  Air  Force  in  image  analysis  operations. 


1.3  Summary  of  Current  Knowledge 

Many  techniques  for  image  restoration  have  been  developed  and  studied 
extensively.  Even  with  these  many  image  restoration  techniques  available,  research 
continues  in  this  area  due  to  limitations  with  each  of  these  techniques.  Many  of  these 
limitations  cause  restored  imagery  to  have  visually  disturbing  effects. 

For  example,  many  techniques  degrade  edges.  Edges  in  an  image  provide 
important  object  definition  but  are  difficult  to  deal  with  in  image  restoration  because  of 
their  high  frequency  content.  The  corruption  of  high  frequencies  causes  blurring  at  these 
edge  locations  in  an  image.  Another  problem  occurring  at  the  location  of  edges  is  the 
Gibbs  phenomenon,  an  interchange  of  overshooting  and  undershooting  a  target  level. 
This  phenomenon  is  associated  with  many  current  signal  and  image  processing 
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techniques.  Thus,  researchers  continually  attempt  to  develop  improved  image  restoration 
techniques. 

One  way  researchers  work  to  solve  these  problems  is  with  adaptive  wavelet 
transforms.  Adaptive  transforms  use  some  or  all  of  the  transform  parameters  to  adapt  to 
the  signal,  typically  to  minimize  some  error  criterion.  The  lifting  construction  [9,  10,  11] 
developed  to  implement  the  wavelet  multiresolution  transform  provides  an  easy  way  to 
introduce  adaptivity  using  prediction  error  energy  as  a  metric.  Adaptable  wavelet 
transforms  have  been  developed  using  the  lifting  construction  [6]  and  show  promise  in 
the  search  for  high  quality  image  restoration. 

A  different  approach  to  image  restoration  was  used  in  [3]  by  using  a  priori 
information  about  an  image.  Perfect  knowledge  of  a  portion  of  the  image  was  assumed 
and  was  used  to  help  remove  noise  from  the  entire  image  in  the  Fourier  transform 
domain.  That  research  showed  that  repeated  use  of  this  information  to  restore  the  image 
using  a  Fourier  transform  technique  improved  the  degraded  image  in  only  the  opposite 
quadrant  of  the  perfect  information.  Other  quadrants  in  the  image  experienced  minimal 
improvement.  To  achieve  these  results,  the  perfect  information  had  to  be  repeatedly 
placed  into  the  image  and  then  the  algorithm  was  run  to  reduce  noise.  The  noise 
reduction  outside  the  area  of  perfect  knowledge  depends  upon  the  noise  correlation 
properties  [3].  For  this  technique,  25%  of  the  image  needs  to  be  perfectly  known  and 
often  knowing  this  much  a  priori  information  about  the  image  is  unrealistic. 
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1.4  Assumptions 

In  this  thesis,  perfect  knowledge  of  an  object  in  a  digital  image  will  be  known. 
Also  the  relative  size  and  rotation  of  this  object  in  the  image  must  be  known.  This 
information  is  necessary  for  development  at  this  stage  due  to  the  fact  that  the  transforms 
used  will  not  be  scale  or  rotation  invariant.  However,  due  to  the  shift  invariant  nature  of 
the  proposed  algorithms,  the  relative  position  of  the  object  can  freely  vary.  Images  will 
be  degraded  with  additive  Gaussian  zero  mean  white  noise.  No  multiplicative  noise  will 
be  introduced  into  the  image  model. 

1.5  Scope 

The  research  will  begin  by  having  complete  knowledge  of  an  object  in  an  image. 
The  orientation,  size  and  intensity  of  the  object  will  be  known.  The  adaptive  wavelet 
transform  will  be  trained  using  this  information;  this  transform  will  then  be  used  for 
image  restoration.  Comparisons  will  be  made  between  these  images  using  soft-threshold 
noise  removal  and  hard-threshold  noise  removal  [2].  Non-adaptive  techniques  will  be 
run  so  that  comparisons  may  be  made.  The  transforms  used  will  be  shift-invariant  so  that 
knowledge  of  the  location  will  not  be  necessary;  this  redundancy  will  also  reduce  the 
pseudo-Gibbs  phenomenon  associated  with  the  non-redundant  discrete  wavelet 
transform-based  restoration  techniques. 
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1.6  Standards 


Images  will  be  evaluated  using  mean  square  error  and  L°°  error  criteria  between 
the  original  image  before  noise  is  added  and  the  image  after  noise  removal.  This 
evaluation  will  be  done  for  the  entire  image  and  also  locally  on  the  object  of  interest. 
More  importantly,  the  visual  appearance  of  the  images  will  be  compared  to  other  image 
restoration  techniques  to  determine  which  is  more  appealing. 


1.7  Approach 

The  multiresolution  transforms  via  lifting  developed  by  Claypoole  [6]  will  be 
modified  to  complete  this  research.  The  transforms  will  be  trained  using  the  known 
information  of  a  local  area  or  object.  Then,  the  trained  transforms  will  be  used  to  process 
an  entire  image  containing  the  training  information  with  degradations.  All  simulations 
will  be  performed  in  MATLAB®. 


II.  Background 


In  its  most  general  form,  image  restoration  is  any  technique  to  recover  an  image 
that  has  been  degraded.  Whether  captured  by  camera,  x-ray,  or  other  optical  systems, 
images  are  typically  degraded  by  various  phenomena.  Movement  of  the  camera  or  other 
external  disruptions  can  cause  blur  in  a  photographic  image,  which  is  visually  disturbing. 
In  digital  imaging,  other  forms  of  interference  degrade  the  image.  In  the  transmission  of 
a  digital  image,  interference  and  thermal  noise  can  cause  degradation  throughout  the 
image.  This  noise  in  an  image  is  purely  an  artifact  of  the  image  process;  one  goal  of 
image  restoration  is  to  remove  that  noise. 

Wavelets  have  become  an  important  tool  for  many  image-processing  techniques. 
They  provide  an  efficient  representation  for  a  broad  range  of  signals.  In  addition  to  this 
parsimonious  representation,  they  provide  spatial  and  frequency  localization,  which  helps 
to  reduce  the  Gibbs  phenomenon  to  only  a  local  effect  hence  referred  to  as  pseudo-Gibbs 
phenomenon  [5].  The  lifting  approach  to  wavelets  provides  an  easy  implementation  and 
a  way  to  introduce  adaptivity.  Redundancy  provides  a  better  representation  by  allowing 
for  any  shift  in  the  data  while  averaging  out  most  of  the  pseudo-Gibbs  phenomenon  and 
other  restoration  artifacts.  This  thesis  will  explore  the  use  of  these  redundant  and 
adaptive  multiresolution  image-processing  techniques  and  a  priori  information  to  enhance 
current  image  restoration  processing  techniques. 
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2.1  Wavelets 


Many  techniques  for  image  restoration  have  been  developed  and  studied 
extensively.  Effective  algorithms  for  removing  degradation  exist  in  the  frequency 
domain.  Mathematical  models  for  image  restoration  are  also  available.  More  recently, 
wavelet-based  image  restoration  techniques  have  proven  to  be  very  effective.  The 
wavelet  decomposition  is  useful  because  signal  information  is  contained  in  a  few 
coefficients.  Another  valuable  feature  of  the  wavelet  decomposition  is  the  fact  that  the 
coefficients  simultaneously  provide  information  on  time/location  and  frequency. 
Wavelet  techniques  preserve  smooth  areas,  polynomials  and  even  piecewise  polynomials 
very  well  but  have  trouble  maintaining  clear,  shaip  edges. 

The  successive  projection  of  a  function  /  into  smaller  orthogonal  subspaces  is  a 
primary  focus  of  wavelet  multiresolution  analysis.  Shifts  and  dilations  of  a  low-pass 
scaling  function  <p(t )  and  a  band  pass  wavelet  function  y/{t)  form  the  smaller  subspaces. 
Presented  below  is  a  derivation  of  the  multiresolution  analysis  for  orthogonal  wavelets. 

Let  {Vmjmsz  be  a  sequence  of  nested  subspaces  in  L2(Z )  and 
Vm  aVm-i  V me  Z  is  required.  Now  let /be  a  function  in  the  subspace  Vm.i  for  some 
m.  The  function  /  e  Vm.j  is  projected  into  the  nested  subspace  Vm  a  Vm.j  by  a 
projection  operator  Pm.  This  projection  operator  Pm  eliminates  the  part  of / which  is  not 
in  Vm  ,  while  the  portion  of  /  in  V„,  remains  undisturbed.  The  orthogonal  projection 
operator,  Qm- 1  -  Pm,  projects  the  function /into  the  subspace  Wm  a  Vm.i.  Here,  I  is  the 
identity  operator.  Qm  and  Pm  are  orthogonal;  therefore,  P,„Qm  -  QmPm  -  0. 
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Let  {JVm]mez  be  another  sequence  of  nested  subspaces  where 
Wm  a  Wm- 1  \/  me  Z .  Let  the  span  of  V„,  u  Wm=  Vm.i  and  let  Vm  and  Wm  be  orthogonal 
subspaces  (i.e.,  Vm-\  =  Vm@Wm  and  VmC\Wm  =  0  ). 

Let  the  set  of  scaling  functions  {0m,i}ie  z  be  an  orthonormal  basis  for  Vm  and  let 

the  set  of  wavelet  functions  {y/m,i\iz  z  be  an  orthonormal  basis  for  Wm.  Vm  c :  Vm.i,  and 
therefore,  0m,  i  can  be  expanded  in  terms  of  0m-i,k  as: 

0m,  /( Jt)  =^Cm-  \{k)0m  -  I,  /t(x).  (2.1) 

k 

Take  the  inner  product  of  0m, i  with  the  basis  function  0m-\,k  to  determine  the 
value  of  the  coefficient  cm.u(k)  as  follows: 

(0m,l,0m-\,k)=  0m,  l(x)0m  -  1,  k{x)dx.  (2.2) 

If  hi,  k  =  ( 0m , /,  0m  - 1,  k) ,  then  the  new  expansion  of  0m, i  becomes: 

0m,l(x )  =  hm  -  \,k0m  -  1,  ^(x).  (2.3) 

k 

Now  to  form  an  orthogonal  wavelet  transform,  the  set  of  functions  0m, i  must 
satisfy  Equation  (2.3)  for  all  m.  For  a  discrete  wavelet  transform,  each  0m, i  must  be 
created  from  integer  shifts  and  dyadic  dilations  (/  and  2m,  respectively)  of  the  scaling 
function  0{x) .  Therefore,  Equation  (2.3)  can  be  written  as: 

0(x-l)  =  '£h(k-2l)0(2x-k).  (2.4) 

k 

Let  g(k  -  21)  =  ( y/m,i,0m-\,k )  and  the  wavelet  equation  can  be  written  as: 

(K*  -/)  =  !>(*-  2X)0{2x  -  k ).  (2.5) 

k 
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Together,  Equations  (2.4)  and  (2.5)  create  the  wavelet  recursion  relations  [8]. 

The  function  f  can  be  completely  described  as  a  linear  combination  of  the  basis 
functions  in  the  subspace  Vm.j.  The  projection  of/into  the  subspace  Vm  can  be  described 
in  terms  of  the  basis  functions  for  Vm.  In  a  similar  way,  the  projection  of  f  into  the 
subspace  Wm  can  be  described  in  terms  of  the  basis  functions  for  Wm.  These  projections 
are: 

[/>»/](*) Cm , 

/ 

[Qnf](x)  =  '£dm,il/fm,l(x), 
l 


where  Cm,  i  =  ( Pmf ,  (pm,  /)  and  dm,  i  =  ( Qmf ,  \j/m,  i) .  The  function / can  be  written  as 

f  —  Pmf  +  Qmf  (2-6) 

because  P+  Q  =  l.  Any  given  decomposition  coefficient  in  V„,  can  be  written  as: 


Cm,  1  —  {Pnf, 

=  {(f-Q4 
~  {f’<Pm'1)- 


(2.7) 


Expand  /  in  terms  of  the  basis  functions  in  Vm.j  and  substitute  the  results  into  Equation 
(2.7)  yielding 


Cm,  /  =  ((£  Cm  ~  - 1  i^J 

—  Cm  -  1,  A  {(pm  -  1,  A,  (pm,  l 


Expressing/in  terms  of  the  basis  functions  for  Wm.j  leads  to: 

dm,  l  —  ^Jcm  -  1,  A  {(pm  -  1,  A, 
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The  inner  products  are  independent  of  the  current  decomposition  level  for  the 
discrete  wavelet  transform.  Therefore, 

=  (2.8) 

(j)m  -\,k)  —  g(k  ~  21).  (2.9) 

The  coefficients  of  the  projection  of  /into  Vm  and  Wm  are  given  by  the  following: 

Cm,,  =  YJCm-\.kh(k-2l),  (2.10) 

k 

rfm,/  =  £cm-Mg(A:-2 /).  (2.11) 

k 

The  coefficients  of /  are  the  starting  point  of  our  decomposition,  and  this  leads  to  a 

filter  bank  implementation  of  the  discrete  wavelet  transform  [8].  As  seen  in  Equation 
2.10  above,  the  coefficients  cmj  are  found  by  convolving  the  coefficients  cm.i,„  with  h, 
and  decimating  the  result  by  2.  The  detail  coefficients  dmJ  are  found  in  similar  fashion 
for  the  filter  g.  These  equations  (2.10  and  2.1 1),  along  with  similar  constructions  for  the 
synthesis  filters,  lead  to  the  filterbank  implementation  of  the  discrete  wavelet  transform, 
shown  in  Figure  2.1. 


x[n] 


Figure  2.1:  Filterbank  Implementation  of  Discrete  Wavelet  Transform 
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After  determining  cmj  and  dmj  ,  the  function  /  can  be  decomposed  onto  the 
orthogonal  subspaces  of  Vm,  i.e.,  V,„+i  and  Wm+i-  The  subspaces  of  Vm  can  be 

decomposed  using  the  same  formulas  because  the  filters  h  and  g  are  independent  of  the 
decomposition  level  m  and  lead  to  a  recursive  decomposition.  The  discrete  wavelet 
transform  consists  of  a  set  of  scaling  coefficients  cm[n],  which  represent  coarse  signal 
information  at  some  lowest  scale  m  =  M,  and  a  set  of  wavelet  coefficients  dm[ri\,  which 
represent  detail  signal  information  at  scales  m  - 

2.2  Denoising 

Denoising  is  an  important  concept  for  image  restoration;  basically,  it  is  any 
method  used  to  remove  noise  from  an  image.  It  has  been  approached  in  many  ways  using 
many  different  techniques.  Denoising  methods  in  the  Fourier  transform  domain  and  in 
the  wavelet  domain  will  be  discussed  but  only  the  latter  will  be  employed  in  this  thesis. 
Although  all  denoising  methods  require  some  a  priori  information  about  the  statistics  of 
the  noise,  denoising  methods  using  a  priori  image  information  will  also  be  reviewed. 

The  process  of  restoring  digital  images  is  often  approached  with  some  a  priori 
knowledge  about  the  cause  of  the  degradation,  typically  the  statistics  of  the  noise.  In 
many  restoration  techniques,  the  degradation  is  modeled  using  this  knowledge  and  then 
the  model  is  used  in  a  way  to  reverse  the  effects  of  the  degradation.  Many  techniques  to 
restore  digital  images  using  this  kind  of  information  have  been  studied  based  on 
frequency  domain  concepts  [1],  These  techniques  typically  involve  some  form  of 
lowpass  or  adaptive  filtering  in  the  frequency  domain.  More  direct  mathematical 
approaches  have  also  been  studied.  Mathematical  approaches  fit  the  data  in  a  least  mean- 
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squares  sense  but  also  require  some  degree  of  smoothness.  In  the  mathematical 
approaches,  the  types  of  a  priori  knowledge  used  in  image  restoration  techniques  vary 
between  equality  and  inequality  constraints.  The  constraints  contain  information  on  the 
cause  of  the  degradation  in  the  image.  For  example,  the  probability  distribution  for  the 
noise  throughout  the  image  could  be  known.  This  image  could  then  be  denoised  using  a 
form  of  low-pass  filter  in  the  Fourier  transform  domain  or  by  using  thresholding  in  the 
wavelet  transform  domain  [2].  The  soft-threshold  in  the  wavelet  domain  provides  a 
smoothness  class  in  the  restored  image  equivalent  to  or  less  than  the  smoothness  class 
held  by  the  original  image  while  maintaining  close  to  a  minimum  in  least  squares  fitting 
to  the  original  image.  Knowledge  about  the  degradation  in  an  image  (in  this  case,  the 
statistics  of  the  noise)  has  been  effectively  used  as  a  priori  information  to  help  restore  the 
image  [2]. 

A  different  approach  to  image  restoration  using  a  priori  information  involves 
specific  knowledge  about  the  image  itself  instead  of  the  degradation.  Using  this 
technique  requires  knowledge  about  the  image  before  the  effects  of  degradation  have 
been  introduced.  A  procedure  in  the  Fourier  transform  domain  using  this  type  of  a  priori 
information  was  developed  in  [3]  and  showed  some  improved  performance.  That 
procedure  requires  full  knowledge  of  one  quarter  of  the  image  without  degradation;  this 
is  not  always  practical.  However,  the  concept  in  general  is  valid;  information  about  the 
image  should  be  a  valuable  tool  for  image  restoration. 

Recently,  the  wavelet  transform  has  shown  great  success  in  reducing  the  noise  in 
images  [2, 4].  A  reconstruction  method  for  an  unknown  function  f  on  [0,1]  with  additive 
noise  was  proposed  with  excellent  results  [2].  The  method  operates  on  the  noisy 
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measured  data  d,  =  /(?,)  +  z,-,  where  i  =  0  1,  /,•  =  //«,  and  z;-  A(0,1).  The  d, 

are  transformed  into  the  wavelet  domain,  and  the  resulting  wavelet  coefficients  are 
translated  towards  zero  by  the  threshold  amount  (2*log(w))1/2*a/(«)1/2  to  form  the 
reconstruction,  ff, .  The  reconstruction,  /*  ,  proved  to  be  as  smooth  as  /  with  high 
probability  in  many  smoothness  measures.  The  estimator  of  the  function  obtained  from 
the  data  fits  the  unknown  function  almost  as  close  in  the  mean  square  sense  as  any 
estimator  could  fit.  These  results  show  how  effective  soft-thresholding  is  for  noise 
removal.  This  is  a  soft-threshold  because  coefficients  with  an  absolute  value  at  or  below 
the  threshold  are  set  to  zero  and  coefficients  with  an  absolute  value  greater  than  the 
threshold  are  reduced  by  the  amount  of  the  threshold.  Hard-thresholding  is  also  available 
in  wavelet  denoising.  The  hard-threshold  is  similar  to  a  soft-threshold  with  the  difference 
being  that  coefficients  with  an  absolute  value  greater  than  the  threshold  are  not  modified. 
However,  hard-thresholding  typically  induces  a  greater  Gibbs  phenomenon  than  soft- 
thresholding  and  the  reconstruction  /*  will  not  be  as  smooth.  Approached  from  a 
mathematics  standpoint,  thresholding  yields  the  best  estimate  [4]. 

Even  with  many  image  restoration  techniques  available,  research  continues  in  this 
area  due  to  problems  with  each  of  these  techniques.  Edges  in  an  image  provide  important 
object  definition  but  are  difficult  to  deal  with  in  image  restoration  because  they  cause 
discontinuities.  Many  image  restoration  techniques  do  not  maintain  the  very  important 
(visually  significant)  edges  of  an  image.  Low  pass  filters  in  the  Fourier  transform 
domain  remove  noise  but  also  cause  sharp  edges  to  blur.  The  Gibbs  phenomenon,  an 
interchange  of  overshooting  and  undershooting  a  target  level,  occurs  globally  in  Fourier- 
based  de-noising.  The  pseudo-Gibbs  phenomenon  appears  only  locally  in  wavelet-based 
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de-noising  [5].  Whether  global  or  local,  the  Gibbs  phenomenon  causes  undesirable  visual 
artifacts  in  the  area  of  discontinuities  after  de-noising  in  the  Fourier  domain  or  the 
wavelet  domain.  Thus,  researchers  continually  attempt  to  improve  image  restoration 
techniques. 

2.3  Redundant  Wavelets  and  Denoising 

Coifman  and  Donoho  developed  a  successful  technique  to  weaken  the  pseudo- 
Gibbs  phenomenon  in  the  wavelet  domain  [5].  “Cycle  spinning”  by  Coifman  “averages 
out”  the  translation  dependence  in  the  wavelet  basis  and  suppresses  the  undesired 
artifacts.  In  this  technique,  the  wavelet  basis  functions  can  be  set  to  shift  a  limited 
number  of  times  or  set  to  shift  over  the  full  scale  of  the  signal,  making  the  transform 
“Fully  Translation-Invariant.”  Either  of  these  techniques  works  well  with  both  hard- 
thresholding  and  soft-thresholding  in  the  wavelet  domain  although  typically  hard 
thresholding  is  used.  In  most  cases,  the  “Fully  Translation-Invariant”  system  gives  the 
best  results  in  terms  of  mean-squared  error  over  “cycle  spinning”  with  a  limited  number 
of  shifts,  although  both  improve  performance  over  non-redundant  denoising.  A 
redundant  undecimated  wavelet  transform  is  also  an  effective  technique  for  improving 
image  restoration.  For  images,  the  wavelet  transform  is  implemented  in  a  separable 
fashion,  i.e.,  the  1-D  wavelet  transform  is  performed  on  the  rows  and  then  the  columns. 
Thus,  for  images  the  redundant  transform  performs  the  wavelet  transform  for  the  image 
with  no  shifts,  one  column  shift,  one  row  shift  and  one  of  each.  This  technique  repeats  at 
each  iteration  so  the  transform  becomes  very  large  quickly  (growing  as  a  factor  of  4L, 
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where  L  is  the  number  of  iterations).  Most  this  information  is  more  than  is  needed  for 
reconstruction  and  can  therefore  be  used  to  improve  restoration. 

2.4  Other  Wavelet-Based  Techniques 

Other  image  processing  techniques  have  been  developed  focusing  on  the 
preservation  of  edges.  In  an  attempt  to  preserve  edges,  Beige  et.  al.  used  wavelet  domain 
image  restoration  and  adaptive  edge-preserving  regularization  [7].  The  distribution  of 
wavelet  coefficients  for  many  images  is  similar  to  that  of  a  Gaussian  density  distribution 
but  the  tails  are  much  heavier  due  to  edges.  This  knowledge  of  the  distribution  for 
wavelet  coefficients  in  an  image  is  used  to  ensure  preservation  of  the  important  edges  in 
images.  The  edge-preserving  image  restoration  techniques  provide  methods  to  remove 
noise  and  degradation  while  maintaining  clear  edges,  but  these  methods  do  have 
problems  in  the  smooth  regions  of  an  image.  This  technique  is  an  example  of  an  adaptive 
denoising  algorithm. 

2.5  Adaptive  Wavelet  Transforms  Via  Lifting 

The  lifting  scheme  [6]  offers  a  simple  implementation  to  perform  the  wavelet 
transform  and  easily  allows  for  perfect  reconstruction  back  to  the  spatial  domain.  The 
lifting  scheme  breaks  the  wavelet  transform  into  a  prediction  and  an  update  step.  Any 
linear,  non-linear,  or  space-varying  function  may  be  used  to  perform  the  predict  and 
update.  The  transform  that  results  is  always  invertible  due  to  the  lifting  construction.  In 
the  linear  prediction  step,  neighbors  of  a  certain  point  are  used  to  predict  that  point.  The 
errors  in  the  prediction  create  the  detail  coefficients  of  the  wavelet  transform.  The 
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prediction  can  be  performed  using  one  point  or  several  points.  These  points  suppress 
polynomials  up  to  the  order  of  the  number  of  points  used  for  prediction.  The  update  step 
uses  the  detail  coefficients  in  the  same  manner  to  update  a  point  and  create  a  coarse 
coefficient  of  the  wavelet  transform.  This  procedure  to  perform  the  wavelet  transform  is 
simple,  logical  and  flexible,  and  will  be  described  in  more  detail  in  Chapter  3. 

Lifting  allows  for  spatially  adaptive  and  scale  adaptive  transforms  [6],  The 
adaptive  transforms  match  the  characteristics  of  a  signal  using  the  data-based  prediction 
measures  of  the  lifting  scheme.  For  the  scale-adaptive  transform  (ScAT),  several  points 
are  used  for  prediction  but  not  all  these  points  are  used  to  suppress  the  polynomials  as  in 
the  usual  lifting  scheme.  These  extra  degrees  of  freedom  are  used  to  adapt  to  the  signal. 
This  adaptation  is  completed  at  each  scale  to  minimize  the  spatially  averaged  squared 
prediction  error.  This  adaptive  predictor  effectively  zeros  out  the  dominant  signal 
structure  at  each  scale.  The  variations  from  this  signal  structure  are  then  represented  by 
the  detail  coefficients  of  the  wavelet  transform.  The  space-adaptive  transform  (SpAT) 
changes  the  wavelet  basis  function  at  each  point  and  scale.  Lifting  allows  the  predictor  to 
change  instantaneously  with  the  signal.  In  this  adaptive  algorithm,  the  predictor  is 
chosen  to  minimize  the  detail  coefficient  of  the  wavelet  transform,  which  is  the  error  in 
the  prediction.  Therefore,  the  number  of  points  used  for  prediction  in  the  space-adaptive 
transform  changes  from  point  to  point  and  scale  to  scale.  The  adaptive  transforms 
perform  almost  as  well  or  even  superior  to  the  common  wavelets  for  denoising.  Again, 
the  available  image  restoration  techniques  are  improved  upon  with  various  algorithms. 
This  thesis  will  concentrate  on  the  ScAT. 
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2.6  Summary 

Research  shows  many  techniques  available  for  image  restoration  and  image  de¬ 
noting.  Certain  factors  in  the  Fourier  and  wavelet  domains  still  cause  problems,  which 
are  visually  undesirable.  Various  methods  improve  some  of  these  effects.  “Cycle 
spinning”  developed  by  Coifman  provides  great  improvement  over  artifacts  caused  by  the 
pseudo-Gibbs  phenomenon  in  the  wavelet  domain  [5].  The  adaptive  transforms 
developed  by  Claypoole  using  the  lifting  technique  can  improve  on  wavelet  denoising  via 
thresholding  [6].  The  use  of  a  priori  knowledge  can  also  improve  the  denoising 
performance.  This  thesis  research  will  combine  wavelets,  adaptivity,  redundancy,  and 
the  use  of  a  priori  information.  The  scale  adaptive  wavelet  transform  will  be  modified  to 
use  redundancy,  then  trained  on  a  known  object  within  an  image.  This  newly  adapted 
transform  will  then  transform  the  entire  degraded  image  into  the  wavelet  domain  where  a 
de-noising  threshold  will  be  applied.  The  inverse  transform  should  yield  improvement  in 
the  area  where  the  known  object  is  located.  A  successful  algorithm  of  this  type  would  be 
desirable  in  many  applications,  particularly  pattern  recognition.  When  looking  for  a 
target  in  pattern  recognition,  an  algorithm  such  as  this  could  be  run  on  an  image  before 
attempting  to  locate  the  target  so  that  the  properties  of  the  target  would  be  enhanced, 
making  the  target  easier  to  identify  by  the  pattern  recognition  system.  This  is  just  one 
example  of  how  this  type  of  locally  adaptive  algorithm  could  be  very  useful  in  image 
restoration  and  image  processing  in  general.  Variations  of  the  adaptive  and  redundant 
wavelet  transforms  have  great  potential  in  enhancing  the  techniques  of  image  restoration. 
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III.  Methodology 


As  described  in  Chapter  II,  Wavelet  transforms  compress  most  signals  into  a  few 
large  coefficients  that  contain  information  on  both  time  and  frequency.  The  lifting 
approach  [9,  10,  11]  to  wavelets  will  be  helpful  in  this  research,  providing  a  structure  to 
easily  incorporate  adaptivity  into  the  wavelet  transform.  The  lifting  interpretation  views 
the  wavelet  transform  as  a  prediction-error  decomposition.  The  lifting  scheme  provides 
easy  implementation  and  a  way  to  introduce  adaptivity.  Redundancy  will  also  be 
implemented  to  reduce  the  psuedo-Gibbs  phenomenon.  Adaptive  transforms  will  be 
utilized  so  that  a  set  of  best  predictors  can  be  found  for  a  given  image;  a  priori 
information  containing  perfect  knowledge  of  part  of  the  image  or  an  object  in  the  image 
will  be  used  to  train  the  adaptive  transform. 

3.1  Lifting 

The  lifting  process  is  reviewed  here.  The  lifting  process  is  performed  in  three 
major  steps:  the  split,  predict,  and  update,  which  are  described  below. 

Split:  The  original  data  is  separated  into  two  subsets.  Any  disjoint  split  of  the 

data  can  be  made  but  the  usual  split  in  the  lifting  scheme  is  the  division  into  even  and  odd 
data  points.  For  the  one-dimensional  data  set  x[«],  the  even  set  would  be  xe[ri\  =  x[2 n] 
and  the  odd  x0[n]  =  x[2n+\}.  Clearly,  the  output  of  the  split  is  time-varying. 
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Predict:  The  error  in  predicting  x0[n]  from  xe[n]  creates  the  wavelet  coefficients 

d[n]  using  the  prediction  operator  P: 

d[n]  =  Xo[n]-P(xe[ri\).  (3.1) 

Update:  The  update  step  forms  a  coarse  approximation  to  the  original  signal  x[n] 

in  the  scaling  coefficients  c[ri\.  They  are  formed  by  combining  xe[n]  and  the  update 

operator  U  applied  to  the  wavelet  coefficients,  d[n\. 

c[n]  =  Xe[n]  +  U(d[n]).  (3.2) 

Multiple  iterations  of  the  lifting  process  are  accomplished  by  repeating  these  steps 
on  the  scaling  coefficients,  c\ri\  to  create  a  complete  set  of  scaling  and  wavelet 
coefficients.  At  each  iteration,  the  scaling  and  wavelet  coefficients  can  be  weighted  to 
normalize  the  energy  of  the  underlying  wavelet  functions. 

Inversion  of  the  lifting  steps  is  easily  accomplished,  regardless  of  the  condition  of 
P  and  U,  by  rearranging  Equations  (3.1)  and  (3.2).  The  forward  and  inverse  lifting 
transforms  use  the  same  P  and  U  to  perfectly  reconstruct  the  transformed  signal.  The 
forward  and  inverse  lifting  stages  are  shown  in  Figures  (3.1)  and  (3.2)  respectively.  Note 
that  these  figures  depict  the  optional  4th  step,  the  weighting  of  the  wavelet  and  scaling 
coefficients  with  ke  and  k0,  respectively. 
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k0 


Figure  3.1:  Forward  Lifting  Stage. 


Figure  3.2:  Inverse  Lifting  Stage. 


In  the  standard  (linear)  lifting  scheme,  an  N-point  predictor  design  is  found  using 
the  desired  predictor’s  z-transform.  For  a  4-point  predictor,  the  z-transform  is 
P(z)  =  pxzA  +P2+  Pi?  +  P4Z2-  An  illustration  to  create  the  wavelet  coefficient  d[ri\  from 
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xe[n],  x0[ti\,  and  P  is  shown  in  Figure  3.3.  By  tracking  the  effect  of  P  up  through  the 
tree,  an  equivalent  prediction  vector  g  can  be  determined.  The  prediction  vector  g  would 
be  applied  to  the  original  data,  x[ri\.  In  this  case,  g  =  [-p\,  0,  -pi,  1,  -pi,  0,  -p^\  . 
Placing  zeros  in  the  locations  corresponding  to  odd  data  points  except  for  the  center  one 
accounts  for  the  split  step  and  allows  the  vector  g  to  be  applied  to  the  original  data. 
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Figure  3.3:  Prediction  Design. 

Now  for  the  N-point  predictor,  an  Nx2N-l  Vandermonde  matrix  V  is  formed. 
The  entries  of  the  Vandermonde  matrix  [V]*,„  =  {n-tf,  where  n  =  1,2,...,2N-1  and  k  = 
0,1,..., N-l  and  /=N  if  a  balanced  filter  is  desired.  Otherwise  /  is  the  number  from  1  to 
2N-1  that  causes  the  column,  where  n-l  =  0,  to  correspond  to  the  1  in  g.  All  even 
columns  are  deleted  to  create  a  new  matrix  V°.  Vg  must  equal  zero  so  that  all  low-order 
polynomials  are  suppressed  by  the  predictor  and  thus  V°p  =  [1  0...0]T.  Now,  p  can  easily 
be  found  by  solving  this  set  of  linear  equations. 
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The  update  filter  uses  an  even  data  point  and  N  of  the  wavelet  coefficients  to 
create  a  coarse  coefficient.  It  is  not  necessary  for  N  to  equal  the  number  of  prediction 
coefficients,  N.  The  prediction  coefficients  are  necessary  to  find  the  update  filter  in  the 
normal  lifting  scheme.  Again,  the  z-transform,  U(z)  =  u\z2  +  ujzA  +  m3  +  114Z,  is  used 
to  find  the  update  filter.  The  update  coefficients  «k  and  the  prediction  coefficients  pk 
combine  to  form  the  update  filter  h  which  should  pass  low-order  polynomials  on  to  the 
coarse  coefficients  c[n ]  and  attenuate  high-order  polynomials.  A  4-point  update 
follows  a  2-point  prediction  in  Figure  3.4  and  tracing  the  effects  of  P  and  U  to  the 
original  data  gives  the  following  update  vector  h: 

T 

h  =  [-pm.  Ml,  -P1U2-P2U1,  U2,  \-p2U2-P\U2,  M3,  -p2U3-p\U4,  M4,  -P2U4]  • 

The  mirror  update  filter  vector  g  where  gn  =  (- 1  )nhn  is  created  to  suppress  low-order 
polynomials. 

g  =  \-p\U\,  -Ml,  -P1U2-P2UU  -u2,  I-P2U2-P1U2,  -M3,  -P2«3-Pl«4,  -U4,  -P2U4]  • 

With  the  previously  determined  prediction  coefficients  p^,  only  the  update  coefficients  Uk 
are  left  to  be  determined.  The  update  coefficients  are  found  by  solving  Vg  =  0,  i.e.,  just 
as  in  the  design  of  the  predictor,  the  mirror  filter  g  is  design  to  suppress  low  order 
polynomials.  Thus,  the  filter  h  will  preserve  these  low  order  polynomials,  yielding  the 
desired  coarse  approximation. 


-22- 


PCUSER1 


PCUSER1 


Microsoft  Word  -  Master  Thesis.doc 
03/09/01  04:18  PM 


x.[0]  x„[0]  x,[1]  x.[1]  x.[2]  x0[2]  x,[3]  x„[3]  x.[4] 


c[2] 


Figure  3.4:  Update  Design. 

With  the  lifting  procedure,  the  transition  to  the  wavelet  domain  is  made 
through  the  steps  of  split,  predict,  and  update.  The  prediction  step  is  designed  to 
eliminate  low-order  polynomial  signal  structures.  The  wavelet  coefficients  formed  by  the 
prediction  step  maintain  the  high-order  polynomials  as  the  details.  Then,  the  update  is 
designed  so  that  the  low-order  polynomial  signal  structure  is  preserved  at  the  next  coarser 
scale. 


3.2  Wavelet  Polyphase  Representation 

The  lifting  wavelet  transform  can  be  written  in  the  form  of  a  perfect 
reconstruction  multirate  filter  bank.  A  perfect  reconstruction  multirate  filter  bank  can  be 
created  using  any  biorthogonal  wavelet  transform. 
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Separating  each  wavelet  filter  into  the  polyphase  components  where 
H(z)  =  He(z2)  +  z"'H0(z2),  the  discrete  wavelet  transform  can  be  completed  as  in  Figure 
3.5. 


Figure  3.5:  Polyphase  Form  of  Wavelet  Transform. 

In  this  form,  down  sampling  the  signal  by  two  and  then  applying  the  filter  He(z)  is  the 
same  as  filtering  using  He(z2)  and  then  down  sampling  by  two  [13]. 

For  the  lifting  scheme  when  the  odd/even  split  is  made,  the  polyphase 
domain  is  created.  In  Figure  3.6,  the  polyphase  matrices  are  created  using  the  predict  and 
update  steps. 
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Figure  3.6:  Lifting  operations  in  Polyphase  Form. 

The  wavelet  coefficients  d[n\  and  the  scaling  coefficients  c[w]  are  created  just  as  before  in 
Equations  (3.1)  and  (3.2).  The  prediction  matrix  does  not  change  the  even  data  points 
x e[ri\  and  the  update  matrix  does  not  change  the  wavelet  coefficients  d[n\.  The  prediction 
matrix  forms  the  wavelet  coefficients  d[n]  as  the  difference  in  the  odd  data  x0[n\  and 
P(xe{n\ ).  The  update  matrix  forms  the  scaling  coefficients  c[n]  as  the  sum  of  even  data 
xe[n\  and  U(d[n]),  the  update  applied  to  the  wavelet  coefficients. 


Combination  of  the  predict  and  update  steps  and  the  undo  update  and  undo  predict 
steps  leads  to  the  lifting  scheme  in  polyphase  form  with  matrices  E(z)  and  R(z)  in  Figure 
3.7. 


Figure  3.7:  Combination  of  Lifting  Steps  in  Polyphase  Form. 
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The  lifting  scheme  typically  does  not  lead  to  orthogonal  transforms.  This  would 

only  be  possible  if  P(z)  and  U(z)  were  constants.  Therefore,  the  lifting  transform  is 

generally  viewed  as  a  biorthogonal  wavelet  transform.  From  the  combination  of  steps  in 

Figure  3.7,  the  wavelet  filters  in  terms  of  the  lifting  scheme  are  as  follows: 

He(z)  =  1 — P(z)U  (z) 

Ho(z)  =  U  (z) 

Ge(Z)  =  -P(z) 

Go(z)  =  1 


Using  the  wavelet  filters  created  from  the  lifting  scheme,  vanishing  moment  constraints 
on  if/{ t)  can  be  mapped  into  constraints  on  the  prediction  and  update  filters  P .  Adding  a 
zero*  vanishing  moment  to  y/(t), 

£^(/)rf/=££g[*M2/-fc)=o 
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and  switching  the  order,  then  £(2r  —  k)dt  is  equal  to  a  constant  mo,  not  zero.  Leaving 
mo£  g [&]={)  and  since  G(z)  =  -P(z2)  +  z~' ,  the  zeroth  vanishing  moment  constraint 

k 

makes 

p\+  p2...+ pN  =  1 . 

Forcing  the  prediction  filter  to  eliminate  zero*  order  polynomials  leads  to  the 
same  constraint.  Additional  vanishing  moments  in  the  analysis  wavelet  function  is 
equivalent  to  eliminating  additional  polynomials  in  the  prediction  step.  After 
determining  the  coefficients  for  the  prediction  filter  or  G(z),  vanishing  moments  are 
added  to  the  synthesis  wavelet  function  vjf(t)  with  the  following 

£  t‘ \j/(t)dt  =  £  X  g[k]tl  $&t-k)  =  0 . 

k 

This  equations  utilizes  the  equivalent  recursion  relation  for  the  dual  (synthesis) 
filters,  as  in  Equation  3.4  above.  Thus,  each  vanishing  moment  on  the  synthesis  wavelet 
function  y/(t)  is  equivalent  to  an  additional  update  filter  polynomial  constraint.  The 
biorthogonal  wavelet  system  created  by  adding  vanishing  moments  to  the  underlying 
wavelet  functions  is  the  same  as  eliminating  and  preserving  polynomials  in  the  predict 
and  update  steps.  Constraints  on  the  wavelet  filters  h,  g,  h,  and  g  are  identical  whichever 
approach  is  used  to  create  them.  The  lifting  scheme  does  not  require  the  use  of  the 
polyphase  matrices  or  the  underlying  scaling  and  wavelet  functions  making  its 
implementation  easier.  This  also  allows  for  nonlinearities  and  adaptivity  to  be  introduced 
into  the  wavelet  transform.  The  properties  of  the  underlying  wavelet  transform  can  be 
controlled  by  forcing  the  prediction  and  update  filters  to  meet  requirements  besides  low- 
order  polynomial  suppression.  In  other  words,  the  transform  can  adapt  by  forcing  the 


-27- 


predict  and  update  filters  to  suppress  the  training  signal,  while  still  maintaining  some 
vanishing  moments  properties  of  the  standard  wavelet  transform. 

3.3  Redundancy 

Superior  performance  in  many  applications  is  gained  by  using  a  redundant 
wavelet  transform  [5].  When  denoising,  the  pseudo-Gibbs  phenomenon  appears  in  the 
local  area  of  discontinuities  after  thresholding  the  non-redundant  wavelet  transform 
coefficients.  The  shift  invariant  or  redundant  transform  improves  thresholding  by 
averaging  over  all  possible  shifts  of  the  signal.  The  redundant  lifting  scheme  is  shown  in 
Figure  3.8  for  a  one-dimensional  signal. 
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This  method  of  shifting  over  all  possible  shifts  and  transforming  is  known  as 
“cycle  spinning”  and  was  developed  by  Coifman  to  “average  out”  translation  dependence 
[5].  Implementation  of  the  redundant  lifted  transform  on  images  will  create  two 
additional  shifts  for  one  iteration,  i.e.,  four  transforms  for  one  image.  As  iterations  L 
increase  the  size  of  the  data  increases  by  4L.  This  causes  a  problem  in  storing  and 
displaying  the  transform.  Instead  of  storing  each  possible  shift  as  a  complete  transform 
(which  would  create  four  transforms  the  same  size  as  the  original  image)  the  transforms 
will  be  stored  as  a  four  dimensional  array.  The  first  two  dimensions  are  the  size  of  the 
original  image,  the  third  is  always  four  for  the  first  four  transforms  created  at  the  first 
iteration,  and  the  fourth  dimension  is  equal  to  the  number  of  iterations.  Thus,  a  transform 
of  a  256x256  image  with  3  iterations  would  have  a  redundant  wavelet  transform 
dimension  of  (256,  256, 4,  3)  and  4  iterations  of  the  same  image  would  have  a  dimension 
of  (256,  256,  4,  4).  For  ease  of  storage  and  display  purposes,  this  method  will  retain  all 
coarse  approximations  at  each  scale,  even  though  only  the  lowest  level  coarse 
approximation  is  needed  for  reconstruction. 

3.4  Adaptivity 

The  lifting  scheme  was  demonstrated  in  Section  3.1  using  polynomial  constraints 
but  restrictions  besides  low-order  polynomial  suppression  can  be  implemented.  This 
research  uses  the  scale-adaptive  transform  to  adapt  the  predictor  of  each  lifting  stage  to 
match  a  training  signal  structure  at  that  scale.  A  N-point  predictor  is  used  to  suppress 
polynomials  up  to  order  M,  less  than  N.  This  is  equivalent  to  adding  M  vanishing 
moments  to  the  underlying  wavelet  function.  The  other  degrees  of  freedom  will  adapt  to 
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the  training  signal.  The  remaining  N-M  degrees  of  freedom  optimize  the  predictor  so 
that  it  minimizes  the  spatially-averaged  squared  prediction  error.  Therefore,  the  predictor 
matches  both  polynomial  and  non-polynomial  signal  structure  at  each  scale.  A  N- 
dimensional  constrained  least-squares  problem  solves  this  optimization  using  the 
constraint  that  the  predictor  must  suppress  M<Nth  order  polynomials.  The  goal  is  to 
minimize  the  sum  of  squared  prediction  errors,  eTe,  when  e  =  Xo  —  X*p.  Here,  x0  are  the 
odd  data  points,  Xe  is  a  matrix  of  the  even  points  such  that  [Xe]n,k=  *e[n-k]  and  p  is  the 
vector  of  prediction  filter  coefficients.  Xe  and  p  are  attempting  to  predict  x0.  Solving 

min||;to-Xp||2  subject  to  V°p  =  [l,0...0]T  results  in  a  p  that  minimizes  eTe,  and 

maintains  M  polynomial  constraints.  In  this  case,  V°  is  an  MxN  matrix  containing  the 
first  M  rows  of  the  full  rank  Vandermonde  matrix  as  determined  in  Section  3.1.  The 
predictor  now  tries  to  capture  the  dominant  training  signal  structure  at  each  scale.  When 
this  transform  is  applied  to  the  degraded  image,  the  noisy  signal  variation  from  the 
training  signal  structure  is  represented  by  the  wavelet  coefficients  </[«].  The  update  filter 
is  designed  so  that  the  coarse  signal  approximation  transformed  at  the  next  scale  should 
maintain  the  low  frequency  dominant  structure  of  the  training  signal. 

3.5  Denoising/Thresholding 

An  adaptive  transform  is  being  implemented  and  thus,  an  appropriate  adaptive 
threshold  will  be  created  for  denoising  so  that  each  sub-band  in  each  iteration  will  have 
its  own  adaptive  threshold.  The  a  priori  information  used  to  train  the  signal  will  also  help 
determine  the  threshold  for  noise  removal.  While  implementing  this  redundant 
transform,  the  transform  must  be  completed  on  the  training  image  as  well  as  the  noisy 
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image  for  each  shift  and  iteration.  This  information  gives  us  not  only  the  training  data  for 
the  next  scale  but  it  also  gives  a  comparison  for  the  size  of  the  wavelet  coefficients  d[n]. 
Looking  separately  at  each  band  a  threshold  is  set  as  the  mean  absolute  value  of  the 
differences  in  the  wavelet  coefficients  for  each  of  the  four  bands  LL,  HL,  LH,  and  HH. 
Thresholds  are  set  for  the  LL  bands  but  these  bands  are  not  thresholded.  The  difference 
in  the  wavelet  coefficients  for  the  training  image  and  the  noisy  image  represent  the  noise 
and  is  a  good  indicator  to  help  establish  an  appropriate  threshold. 
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IV.  Results 


In  this  chapter,  the  success  of  denoising  with  the  redundant  adaptive  wavelet  transform 
using  a  priori  information  is  demonstrated.  As  discussed  in  Chapter  2,  current  image 
restoration  techniques  have  certain  undesirable  effects.  In  this  research,  a  new  technique, 
combining  adaptive  wavelets,  redundancy,  and  a  priori  information  was  developed.  This 
thesis  research  gives  reasons  why  such  an  algorithm  should  be  powerful  to  image- 
processing  techniques.  Now,  the  performance  of  this  algorithm  is  compared  to  other 
redundant  wavelet-based  techniques. 


4.1  Error  Measures 

Two  errors  measures  are  used  to  measure  the  performance  of  the  adaptive 
algorithm.  The  peak  signal-to-noise  ratio  (PSNR)  is  used  to  help  determine  the  strength 
of  the  true  signal  compared  to  the  reconstructed  signal  (or  noisy  signal). 

f  ^ 


PSNR  =20*  log  10 


max  Xi 


V  N 


(4.1) 


where  x ,■  is  the  Ith  pixel  of  the  original  image,  x,  is  the  /th  pixel  of  the  reconstructed  (or 
noisy)  image,  and  N  is  the  total  number  of  pixels.  The  PSNR  values  are  determined 
before  noise  removal  and  after  for  comparison  of  the  improvement  made  in  the 
reconstructed  image.  PSNR  is  a  common  metric  in  the  image  processing  community. 
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The  second  measure  of  performance  used  for  the  algorithm  is  L“.  The  error  here 
is  the  maximum  difference  in  the  desired  image  and  the  denoised  image.  This  metric 
relates  to  the  human  visual  system,  and  is  a  good  indicator  of  the  visual  “quality”  of  an 
image.  These  performance  measures  are  good  indicators  of  how  the  training  algorithm 
compares  to  other  types  of  wavelet  denoising. 

4.2  Denoising  Experiments 

Several  wavelet  transforms  are  used  to  transform  an  image.  The  image  is  then 
denoised  while  in  the  wavelet  domain  before  it  is  reconstructed.  For  a  basis  comparison, 
a  (4,4)  lifting  scheme  is  used  with  4  vanishing  moments  so  that  the  transform  is  not 
adaptive  and  does  not  take  advantage  of  any  a  priori  information.  Also  for  comparison,  a 
separate  baseline  case  is  completed.  In  this  case,  the  (4,4)  lifting  transform  is  used  with 
only  2  vanishing  moments  and  this  transform  adapts  to  the  noisy  image.  In  other  words, 
this  transform  is  not  trained  using  the  a  priori  information  but  it  does  use  adaptivity.  For 
each  of  these  transforms,  a  soft  threshold  was  used  with  the  threshold  being  3  times  the 
standard  deviation  of  the  noise.  A  hard  threshold  was  also  completed  on  these  transforms 
using  3.6  times  the  standard  deviation  of  the  noise  as  the  threshold.  These  are  the 
standard  thresholds  experimentally  derived  from  [12]. 

One  test  of  the  training  algorithm  is  completed  using  a  (4,4)  lift  with  2  vanishing 
moments  adapted  to  the  training  image  and  for  consistency  to  the  baseline  comparisons,  a 
soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise  and  a  hard  threshold 
equal  to  3.6  times  the  standard  deviation  of  the  noise.  A  second  test  of  the  algorithm  also 
uses  a  (4,4)  lift  with  2  vanishing  moments  adapted  to  the  training  image  and  for  denoising 
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uses  the  threshold  set  adaptively  by  the  training  algorithm  as  discussed  in  Chapter  3. 
This  threshold  was  multiplied  by  various  constants  to  find  which  produced  the  best 
results  and  then  used  for  denoising.  With  a  soft  threshold,  the  adaptive  threshold  is 
multiplied  by  a  constant  of  1.5  for  the  best  results.  When  a  hard  threshold  is  used,  the 
adaptive  threshold  is  multiplied  by  a  constant  of  3  for  the  best  results.  So  one  test  uses 
the  a  priori  information  to  train  the  transform  but  uses  a  standard  threshold  and  the 
second  test  uses  the  a  priori  information  to  train  the  transform  and  to  set  the  adaptive 
threshold. 

The  original  image  is  made  zero  mean  and  is  show  in  Figure  4.1.  To  obtain  the 
training  image  in  Figure  4.2,  a  part  of  the  original  image  is  made  zero  mean  and  then 
padded  with  zeros  so  that  it  is  the  same  size  as  the  test  image.  This  is  necessary  because 
the  algorithm  calls  for  an  equal  number  of  iterations  of  the  wavelet  transform  to  be 
completed  on  both  the  original  image  and  the  training  image  and  therefore  to  make  this 
possible  at  all  times,  the  images  must  be  the  same  size. 

4.3  Denoising  Results 

The  first  noisy  image  in  Figure  4.3  is  created  by  adding  a  noise  amount  of 
approximately  25  times  the  standard  deviation  of  the  noise  giving  a  PSNR  value  of 
approximately  20  dB  and  a  L°°  measurement  of  105.38.  This  represents  a  high  amount  of 
noise  for  the  image. 

Results  after  using  a  soft  threshold  are  looked  at  first.  The  image  reconstructed 
after  using  a  non-adaptive  transform  and  a  standard  threshold  has  significant  blurring 
around  the  edges  of  the  image,  but  the  PSNR  is  improved  by  approximately  8  dB.  This 
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first  baseline  case  is  shown  in  Figure  4.4.  The  second  baseline,  in  Figure  4.5,  case 
improves  the  PSNR  and  the  L°°  measures  for  the  whole  image  but  blurring  is  still  present. 

Again  looking  at  the  soft  threshold  results,  the  test  using  an  adaptive  threshold 
and  the  training  algorithm  reduces  blurring  but  effects  of  the  noise  are  still  present,  see 
Figure  4.6.  The  PSNR  values  are  comparable  to  the  baseline  measurements  and  the  L°° 
values  are  reduced.  Using  the  training  algorithm  with  the  standard  threshold  (Figure  4.7) 
produces  an  image  that  is  visually  as  good  as  the  first  baseline  case,  possibly  even  a  little 
less  blur.  The  error  values  are  not  significantly  reduced  but  are  comparable.  All  soft 
threshold  results  with  this  high  noise  case  cause  some  blurring  of  the  edges. 

Most  hard  thresholds  improve  over  their  equivalent  soft  thresholds.  The  hard 
threshold  results  are  visibly  better  because  there  is  less  blurring.  In  Figure  4.8,  the  first 
baseline  case  is  an  improvement  over  previous  transforms  but  the  edges  still  acquire 
undesired  effects.  In  Figure  4.9,  adapting  to  the  noisy  signal  leaves  obvious  artifacts  in 
the  reconstructed  image  but  edge  effects  are  not  as  obvious.  Artifacts  left  from  the  noise 
are  not  as  prominent  in  Figure  4.10.  This  reconstructed  image  has  less  blur,  yet  it  is  also 
smoother.  The  best  PSNR  value  for  the  high  noise  case  is  produced  by  Figure  4.1 1  where 
the  training  algorithm  is  used  with  the  standard  threshold.  This  test  also  gives  a  good  L°° 
value. 

The  training  error  measurements  were  very  similar  to  the  overall  picture 
measurements.  This  was  probably  because  the  majority  of  the  edges  in  the  image  where 
within  the  training  area  without  any  large  smooth  areas  as  those  found  outside  of  this 
training  area. 
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Figure  4. 1 :  Original  star  image  before  the  addition  of  noise. 
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Figure  4.2:  Training  star  image  used  for  adapting  the  transform. 


Figure  4.3:  Star  image  and  additive  Gaussian  noise  with  standard  deviation  equal  to 
10%  of  the  maximum  signal  value  (high  noise  case).  For  entire  image,  PSNR  =  20.15 
dB  and  L“  =  105.38.  In  the  training  area,  PSNR  =  20.12  dB  and  L°°  =  105.38. 
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Figure  4.4:  Star  image  denoised  using  (4,4)  lift,  4  vanishing  moments  (non-adaptive) 

and  a  soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise.  For  entire 
image,  PSNR  =  28.00  dB  and  L°°  =  84.96.  In  the  training  area,  PSNR  =  25.86  dB  and 
L“  =66.14. 


Figure  4.5:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
noisy  signal  and  a  soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise.  For 
entire  image,  PSNR  =  29.09  dB  and  L“  =  80.16.  In  the  training  area, 
PSNR  =  27.15  dB  and  L~  =  67.20. 


Figure  4.6:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
training  signal  and  a  soft  threshold  equal  to  1.5  times  the  adaptive  threshold  set  by  the 
transform.  For  entire  image,  PSNR  =  27.44  dB  and  L°°  =  76.67.  In  the  training  area, 
PSNR  =  26.03  dB  and  L~  =  66.02. 
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Figure  4.7:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
training  signal  and  a  soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise. 
For  entire  image,  PSNR  =  28.43  dB  and  L”  =  86.09.  In  the  training  area, 
PSNR  =  26.33  dB  and  L“  =  68.66. 


Figure  4.8:  Star  image  denoised  using  (4,4)  lift,  4  vanishing  moments  (non-adaptive) 

and  a  hard  threshold  equal  to  3.6  times  the  standard  deviation  of  the  noise.  For  entire 
image,  PSNR  =  29.91  dB  and  L”  =  71.68.  In  the  training  area,  PSNR  =  27.73  dB  and 
L“  =64.80. 


Figure  4.9:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
noisy  signal  and  a  hard  threshold  equal  to  3.6  times  the  standard  deviation  of  the  noise. 
For  entire  image,  PSNR  =  28.82  dB  and  V°  =  75.20.  In  the  training  area, 
PSNR  =  27.39  dB  and  L“  =  75.20. 


The  second  noisy  image  in  Figure  4.12  is  created  by  adding  a  noise  amount  of 
approximately  15  times  the  standard  deviation  of  the  noise  giving  a  PSNR  value  of 
approximately  25  dB  and  a  L“  measurement  of  105.38.  This  represents  a  medium 
amount  of  noise  for  the  image. 

Soft  threshold  results  are  looked  at  first.  The  image  reconstructed  after  using  a 
non-adaptive  transform  and  a  standard  threshold  has  blurring  around  the  edges  of  the 
image,  but  the  PSNR  is  improved  by  approximately  5  dB.  This  first  baseline  case  is 
shown  in  Figure  4.13.  The  second  baseline,  in  Figure  4.14,  case  improves  the  PSNR  and 
the  L“  measures  for  the  whole  image  but  blurring  is  still  present. 

Again  looking  at  the  soft  threshold  results,  the  test  using  an  adaptive  threshold 
and  the  training  algorithm  reduces  blurring  but  effects  of  the  noise  are  still  present,  see 
Figure  4.15.  The  PSNR  values  are  comparable  to  the  baseline  measurements  and  the  L“ 
values  are  significantly  reduced.  Using  the  training  algorithm  with  the  standard  threshold 
(Figure  4.16)  produces  an  image  that  is  visually  as  good  as  the  first  baseline  case,  and  has 
some  improvements  in  the  measurements. 

Most  hard  thresholds  improve  over  their  equivalent  soft  thresholds.  The  hard 
threshold  results  are  visibly  better  because  there  is  less  blurring.  In  Figure  4.17,  the  first 
baseline  case  is  an  improvement  over  previous  transforms  but  the  edges  still  blur.  In 
Figure  4.18,  adapting  to  the  noisy  signal  leaves  obvious  artifacts  in  the  reconstructed 
image  but  edges  are  clearer.  Artifacts  left  from  the  noise  still  exist  but  are  reduced  in  size 
for  Figure  4.19.  This  adaptive  hard  threshold  used  with  the  training  algorithm  gives 
excellent  L°°  values  in  this  example.  Figure  4.20  produces  the  best  PSNR  value  for  the 
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medium  noise  case  which  is  where  the  training  algorithm  is  used  with  the  standard 
threshold.  This  test  also  gives  a  good  L°°  value. 

The  training  error  measurements  were  very  similar  to  the  overall  picture 
measurements.  This  was  probably  because  the  majority  of  the  edges  in  the  image  where 
within  the  training  area  without  any  large  smooth  areas  as  those  found  outside  of  this 
training  area. 


-50- 


Figure  4.12:  Star  image  and  additive  Gaussian  noise  with  standard  deviation  equal  to 
6%  of  the  maximum  signal  value  (medium  noise  case).  For  entire  image,  PSNR  =  24.58 
dB  and  L°°  =  63.23.  In  the  training  area,  PSNR  =  24.56  dB  and  L“  =  63.23. 


Figure  4.13:  Star  image  denoised  using  (4,4)  lift,  4  vanishing  moments  (non-adaptive) 
and  a  soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise.  For  entire 
image,  PSNR  =  29.88  dB  and  L“  =  70.87.  In  the  training  area,  PSNR  =  27.69  dB  and 
L“  =58.83.  (medium  noise) 


Figure  4.14:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
noisy  signal  and  a  soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise.  For 
entire  image,  PSNR  =  30.91  dB  and  L“  =  66.62.  In  the  training  area, 
PSNR  =  28.81  dB  and  L~  =  55.76. 


Figure  4.15:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
training  signal  and  a  soft  threshold  equal  to  1 .5  times  the  adaptive  threshold  set  by  the 
transform.  For  entire  image,  PSNR  =  29.92  dB  and  L“  =  57.74.  In  the  training  area, 
PSNR  =  28.34  dB  and  IT  =  50.86. 
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Figure  4. 16:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
training  signal  and  a  soft  threshold  equal  to  3  times  the  standard  deviation  of  the  noise. 
For  entire  image,  PSNR  =  30.27  dB  and  L°°  =  72.45.  In  the  training  area, 

PSNR  =  28.12  dB  and  L“  =  59.20. 


Figure  4.17:  Star  image  denoised  using  (4,4)  lift,  4  vanishing  moments  (non-adaptive) 
and  a  hard  threshold  equal  to  3.6  times  the  standard  deviation  of  the  noise.  For  entire 
image,  PSNR  =  32.32  dB  and  L°°  =  54.04.  In  the  training  area,  PSNR  =  30.22  dB  and 
L°°  =54.04. 
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Figure  4. 18:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
noisy  signal  and  a  hard  threshold  equal  to  3.6  times  the  standard  deviation  of  the  noise. 
For  entire  image,  PSNR  =  31.93  dB  and  L“  =  56.08.  In  the  training  area, 
PSNR  =  30.18  dB  and  L“  =  56.08. 


Figure  4.19:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
training  signal  and  a  hard  threshold  equal  to  3  times  the  adaptive  threshold  set  by  the 
transform.  For  entire  image,  PSNR  =  31.59  dB  and  L°°  =52.10.  In  the  training  area, 
PSNR  =  30.19  dB  and  L“  =51.81. 


Figure  4.20:  Star  image  denoised  using  (4,4)  lift,  2  vanishing  moments,  adaptive  to 
training  signal  and  a  hard  threshold  equal  to  3.6  times  the  standard  deviation  of  the  noise. 
For  entire  image,  PSNR  =  32.33  dB  and  L“  =  57.43.  In  the  training  area, 
PSNR  =  30.28  dB  and  L“  =  57.43. 
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Table  4.3 :  PSNR  measurements  for  medium  noise  case. 
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4.4  Summary 

The  new  algorithm  used  with  the  adaptive  soft  threshold  performs  very  well  by  L“ 
measures  and  comparable  by  PSNR  measures  for  both  the  whole  image  and  the  training 
area.  The  new  algorithm  used  with  the  non-adaptive  hard  threshold  gives  the  best  overall 
PSNR  measure,  again  for  both  the  whole  image  and  the  training  area.  In  each 
thresholding  case,  the  new  algorithm  is  better  than  or  comparable  to  the  (4,4)  lift  with  4 
vanishing  moments  and  the  (4,4)  lift  with  2  vanishing  moments  adapted  to  the  noisy 
signal  in  the  PSNR  and  L“  measures. 

Visually,  the  training  algorithm  results  are  very  good.  Images  resulting  from  the 
training  algorithms  using  soft-thresholds  have  fewer  artifacts  than  the  other  images  after  a 
soft  threshold  was  used.  The  adaptive  threshold  used  with  hard  thresholding  did  not 
perform  as  well  in  error  measures  but  the  image  is  smoother  with  fewer  artifacts  than  the 
other  images  after  hard  thresholding. 

These  results  show  that  the  training  algorithm  performs  very  well  in  the  area  of 
training  and  sometimes  the  entire  image.  The  training  algorithm  obtains  the  best  PSNR 
value  overall  the  images  for  the  whole  image  and  the  training  area.  The  combination  of 
the  adaptive  wavelets,  redundancy  and  a  priori  information  have  proven  to  be  useful  in 
obtaining  better  image  restoration. 
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V.  Conclusion 


This  research  explored  many  techniques  available  for  image  restoration  and  image 
denoising.  Certain  factors  in  the  Fourier  and  wavelet  domains  still  cause  problems, 
which  are  visually  undesirable.  As  technology  progresses,  more  improvement  in  image 
restoration  is  desired,  and  can  be  achieved  by  incorporating  redundancy  and  adaptivity 
into  the  wavelet  transform.  In  this  thesis  research,  a  Scale  Adaptive  Transform  (ScAT) 
was  extended  to  form  a  redundant  transform  for  a  2-dimensional  signal.  The  ScAT  was 
further  modified  so  that  it  adapts  to  a  training  signal  and  then  transforms  a  separate 
signal.  In  this  research,  the  modified  ScAT  is  trained  on  an  image  containing  a  priori 
information,  typically  a  known  object  in  the  image.  The  noisy  image  is  then  transformed 
into  the  wavelet  domain  using  the  newly  trained  transform,  where  various  denoising 
threshold  techniques  were  applied  to  remove  noise  degradations.  These  techniques 
improved  the  image  in  the  area  of  the  known  object  and  at  times  even  improved  the 
whole  image.  The  algorithms  and  results  presented  in  this  thesis  demonstrate  that 
improved  denoising  is  possible  with  the  incorporation  of  redundancy  and  adaptivity  into 
existing  wavelet-based  image  restoration  techniques. 

5.1  Applications 

An  algorithm  of  this  type  could  be  used  as  a  pre-processing  step  for  many  image- 
processing  applications.  In  target  recognition,  researchers  look  for  distinguishing 
features  of  a  possible  target  in  an  image  to  allow  a  computer  to  locate  and  identify  the 
target.  An  algorithm  that  enhances  the  features  of  the  target  would  greatly  improve  the 
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probability  of  its  detection  in  a  cluttered  image.  This  is  just  one  example  of  how  this 
type  of  redundant  adaptive  algorithm  could  be  very  useful  in  image  restoration  and  image 
processing.  In  general,  any  image  analysis  application  where  a  priori  information  is 
known  about  the  image  could  be  improved  by  the  algorithms  presented  in  this  thesis. 

5.2  Recommendations  for  Future  Research 

The  success  of  this  algorithm  gives  encouragement  to  continue  this  area  of 
research.  This  technique  improves  noise  removal  in  an  image  using  redundancy,  a  priori 
information  and  just  one  adaptive  prediction  filter  per  scale;  therefore,  if  the  prediction 
filters  also  adapted  spatially  (e.g.,  multiple  filters  per  scale),  the  improvement  could  be 
even  greater.  This  modification  would  require  developing  a  method  to  compensate  with 
appropriate  filters  outside  of  the  training  area  where  only  zeros  are  located  (training  on 
data  comprised  of  all  zeros  is  an  ill-posed  problem).  Also,  additional  storage  would  be 
required  for  these  additional  prediction  and  update  filters.  This  is  a  difficult  task,  due  to 
the  nature  of  the  redundant  wavelet  transform.  However,  the  creation  of  such  a  redundant 
spatially  adaptive  transform  would  obviously  further  enhance  image  restoration 
techniques. 
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